home *** CD-ROM | disk | FTP | other *** search
/ STraTOS 1997 April & May / STraTOS 1 - 1997 April & May.iso / CD01 / INTERNET / SITES / LITTLE / P3SRC.ZIP / ATARI / LATHE.C < prev    next >
Encoding:
C/C++ Source or Header  |  1996-04-25  |  31.4 KB  |  1,393 lines
  1. /****************************************************************************
  2. *                   lathe.c
  3. *
  4. *  This module implements functions that manipulate lathes.
  5. *
  6. *  This module was written by Dieter Bayer [DB].
  7. *
  8. *  from Persistence of Vision(tm) Ray Tracer
  9. *  Copyright 1996 Persistence of Vision Team
  10. *---------------------------------------------------------------------------
  11. *  NOTICE: This source code file is provided so that users may experiment
  12. *  with enhancements to POV-Ray and to port the software to platforms other
  13. *  than those supported by the POV-Ray Team.  There are strict rules under
  14. *  which you are permitted to use this file.  The rules are in the file
  15. *  named POVLEGAL.DOC which should be distributed with this file. If
  16. *  POVLEGAL.DOC is not available or for more info please contact the POV-Ray
  17. *  Team Coordinator by leaving a message in CompuServe's Graphics Developer's
  18. *  Forum.  The latest version of POV-Ray may be found there as well.
  19. *
  20. * This program is based on the popular DKB raytracer version 2.12.
  21. * DKBTrace was originally written by David K. Buck.
  22. * DKBTrace Ver 2.0-2.12 were written by David K. Buck & Aaron A. Collins.
  23. *
  24. *****************************************************************************/
  25.  
  26. /****************************************************************************
  27. *
  28. *  Explanation:
  29. *
  30. *    The lathe primitive is defined by a set of points in 2d space which
  31. *    are interpolated by linear, quadratic, or cubic splines. The resulting
  32. *    2d curve is rotated about an axis to form the final lathe object.
  33. *
  34. *    All calculations are done in the object's (u,v,w)-coordinate system
  35. *    with the (w)-axis being the rotation axis.
  36. *
  37. *    One spline segment in the (r,w)-plane is given by the equations
  38. *
  39. *      fw(t) = Aw * t^3 + Bw * t^2 + Cw * t + Dw  and
  40. *      fr(t) = Ar * t^3 + Br * t^2 + Cr * t + Dr,
  41. *
  42. *    with the parameter t ranging from 0 to 1 and r = sqrt(u*u+v*v).
  43. *
  44. *    To intersect a ray R = P + k * D transformed into the object's
  45. *    coordinate system with the lathe object, the equations
  46. *
  47. *      (Pu + k * Du)^2 + (Pv + k * Dv)^2 = (fr(t))^2
  48. *                          (Pw + k * Dw) = fw(t)
  49. *
  50. *    have to be solved for t. For valid intersections (0 <= t <= 1)
  51. *    the corresponding k can be calculated from one of the above equations.
  52. *
  53. *    Note that the degree of the polynomal to solve is two times the
  54. *    degree of the spline segments used.
  55. *
  56. *    Note that Pu, Pv, Pw and Du, Dv, Dw denote the coordinates
  57. *    of the vectors P and D.
  58. *
  59. *  Syntax:
  60. *
  61. *    lathe
  62. *    {
  63. *      [ linear_spline | quadratic_spline | cubic_spline ]
  64. *
  65. *      number_of_points,
  66. *
  67. *      <P[0]>, <P[1], ..., <P[n-1]>
  68. *
  69. *      [ sturm ]
  70. *    }
  71. *
  72. *    Note that the P[i] are 2d vectors.
  73. *
  74. *    Linear interpolation is used by default. In this case all n Points
  75. *    are used. In the quadratic case the first point is used to determine
  76. *    the derivates at the starting point P[1]. In the cubic case
  77. *    the first and last points are used to determine the derivatives at
  78. *    the starting point P[1] and ending point P[n-2].
  79. *
  80. *    To get a closed (and smooth) curve you have make sure that
  81. *
  82. *      P[0] = P[n-1] in the linear case,
  83. *
  84. *      P[0] = P[n-2] and P[1] = P[n-1] in the quadratic case, and
  85. *
  86. *      P[0] = P[n-3] and P[1] = P[n-2] and P[2] = P[n-1] in the cubic case.
  87. *
  88. *    Note that the x coordinate of a point corresponds to the r coordinate
  89. *    and the y coordinate to the w coordinate;
  90. *
  91. *  ---
  92. *
  93. *  Ideas for the lathe were taken from:
  94. *
  95. *    P. Burger and D. Gillies, "Rapid Ray Tracing of General Surfaces
  96. *    of Revolution", New Advances in Computer Graphics, Proceedings
  97. *    of CG International '89, R. A. Earnshaw, B. Wyvill (Eds.),
  98. *    Springer, ..., pp. 523-531
  99. *
  100. *    P. Burger and D. Gillies, "Swept Surfaces", Interactive Computer
  101. *    Graphics, Addison-Wesley, 1989, pp. 376-380
  102. *
  103. *  ---
  104. *
  105. *  Jun 1994 : Creation. [DB]
  106. *
  107. *****************************************************************************/
  108.  
  109. #include "frame.h"
  110. #include "povray.h"
  111. #include "vector.h"
  112. #include "povproto.h"
  113. #include "bbox.h"
  114. #include "lathe.h"
  115. #include "polysolv.h"
  116. #include "matrices.h"
  117. #include "objects.h"
  118. #include "torus.h"
  119.  
  120.  
  121.  
  122. /*****************************************************************************
  123. * Local preprocessor defines
  124. ******************************************************************************/
  125.  
  126. /* Minimal intersection depth for a valid intersection. */
  127.  
  128. #define DEPTH_TOLERANCE 1.0e-4
  129.  
  130. /* |Coefficients| < COEFF_LIMIT are regarded to be 0. */
  131.  
  132. #define COEFF_LIMIT 1.0e-20
  133.  
  134. /* Max. number of intersecions per spline segment. */
  135.  
  136. #define MAX_INTERSECTIONS_PER_SEGMENT 6
  137.  
  138.  
  139.  
  140. /*****************************************************************************
  141. * Local typedefs
  142. ******************************************************************************/
  143.  
  144. typedef struct Lathe_Intersection_Structure LATHE_INT;
  145.  
  146. struct Lathe_Intersection_Structure
  147. {
  148.   DBL d;  /* Distance of intersection point                  */
  149.   DBL w;  /* Paramter of intersection point with n-th spline */
  150.   int n;  /* Number of segment hit                           */
  151. };
  152.  
  153.  
  154.  
  155. /*****************************************************************************
  156. * Static functions
  157. ******************************************************************************/
  158.  
  159. static int intersect_lathe PARAMS((RAY *Ray, LATHE *Lathe, LATHE_INT *Int));
  160. static int   All_Lathe_Intersections PARAMS((OBJECT *Object, RAY *Ray, ISTACK *Depth_Stack));
  161. static int   Inside_Lathe PARAMS((VECTOR point, OBJECT *Object));
  162. static void  Lathe_Normal PARAMS((VECTOR Result, OBJECT *Object, INTERSECTION *Inter));
  163. static void  *Copy_Lathe PARAMS((OBJECT *Object));
  164. static void  Translate_Lathe PARAMS((OBJECT *Object, VECTOR Vector, TRANSFORM *Trans));
  165. static void  Rotate_Lathe PARAMS((OBJECT *Object, VECTOR Vector, TRANSFORM *Trans));
  166. static void  Scale_Lathe PARAMS((OBJECT *Object, VECTOR Vector, TRANSFORM *Trans));
  167. static void  Transform_Lathe PARAMS((OBJECT *Object, TRANSFORM *Trans));
  168. static void  Invert_Lathe PARAMS((OBJECT *Object));
  169. static void  Destroy_Lathe PARAMS((OBJECT *Object));
  170.  
  171.  
  172. /*****************************************************************************
  173. * Local variables
  174. ******************************************************************************/
  175.  
  176. static METHODS Lathe_Methods =
  177. {
  178.   All_Lathe_Intersections,
  179.   Inside_Lathe, Lathe_Normal,
  180.   Copy_Lathe,
  181.   Translate_Lathe, Rotate_Lathe,
  182.   Scale_Lathe, Transform_Lathe, Invert_Lathe, Destroy_Lathe
  183. };
  184.  
  185.  
  186.  
  187. /*****************************************************************************
  188. *
  189. * FUNCTION
  190. *
  191. *   All_Lathe_Intersections
  192. *
  193. * INPUT
  194. *
  195. *   Object      - Object
  196. *   Ray         - Ray
  197. *   Depth_Stack - Intersection stack
  198. *   
  199. * OUTPUT
  200. *
  201. *   Depth_Stack
  202. *   
  203. * RETURNS
  204. *
  205. *   int - TRUE, if a intersection was found
  206. *   
  207. * AUTHOR
  208. *
  209. *   Dieter Bayer
  210. *   
  211. * DESCRIPTION
  212. *
  213. *   Determine ray/lathe intersection and clip intersection found.
  214. *
  215. * CHANGES
  216. *
  217. *   Jun 1994 : Creation.
  218. *
  219. ******************************************************************************/
  220.  
  221. static int All_Lathe_Intersections(Object, Ray, Depth_Stack)
  222. OBJECT *Object;
  223. RAY *Ray;
  224. ISTACK *Depth_Stack;
  225. {
  226.   int i, imax, Found;
  227.   LATHE_INT *I;
  228.   VECTOR IPoint;
  229.  
  230.   /* Allocate memory for intersections (two per segment). */
  231.  
  232.   I = (LATHE_INT *)POV_MALLOC(MAX_INTERSECTIONS_PER_SEGMENT*((LATHE *)Object)->Number*sizeof(LATHE_INT), "lathe intersection array");
  233.  
  234.   Found = FALSE;
  235.  
  236.   imax = intersect_lathe(Ray, (LATHE *)Object, I);
  237.  
  238.   if (imax)
  239.   {
  240.     for (i = 0; i < imax; i++)
  241.     {
  242.       if ((I[i].d > DEPTH_TOLERANCE) && (I[i].d < Max_Distance))
  243.       {
  244.         VEvaluateRay(IPoint, Ray->Initial, I[i].d, Ray->Direction);
  245.  
  246.         if (Point_In_Clip(IPoint, Object->Clip))
  247.         {
  248.           push_entry_i1_d1(I[i].d,IPoint,Object,I[i].n,I[i].w,Depth_Stack);
  249.  
  250.           Found = TRUE;
  251.         }
  252.       }
  253.     }
  254.   }
  255.  
  256.   POV_FREE (I);
  257.  
  258.   return(Found);
  259. }
  260.  
  261.  
  262.  
  263. /*****************************************************************************
  264. *
  265. * FUNCTION
  266. *
  267. *   intersect_lathe
  268. *
  269. * INPUT
  270. *
  271. *   Ray          - Ray
  272. *   Lathe        - Lathe
  273. *   Intersection - Lathe intersection structure
  274. *   
  275. * OUTPUT
  276. *
  277. *   Intersection
  278. *   
  279. * RETURNS
  280. *
  281. *   int - Number of intersections found
  282. *   
  283. * AUTHOR
  284. *
  285. *   Dieter Bayer
  286. *   
  287. * DESCRIPTION
  288. *
  289. *   Determine ray/lathe intersection.
  290. *
  291. *   NOTE: The curve is rotated about the y-axis!
  292. *         Order reduction cannot be used.
  293. *
  294. * CHANGES
  295. *
  296. *   Jun 1994 : Creation.
  297. *
  298. ******************************************************************************/
  299.  
  300. static int intersect_lathe(Ray, Lathe, Intersection)
  301. RAY *Ray;
  302. LATHE *Lathe;
  303. LATHE_INT *Intersection;
  304. {
  305.   int i, ii, j, n1, n2, skip;
  306.   DBL k, len, r, m, w, Dy2, r0;
  307.   DBL x1[7], x2[3], y1[6], y2[2];
  308.   VECTOR P, D;
  309.   LATHE_SPLINE_ENTRY Entry;
  310.  
  311.   Increase_Counter(stats[Ray_Lathe_Tests]);
  312.  
  313.   /* Init intersection structure. */
  314.  
  315.   Intersection->d = 
  316.   Intersection->w = 0.0;
  317.   Intersection->n = 0;
  318.  
  319.   /* Transform the ray into the lathe space. */
  320.  
  321.   MInvTransPoint(P, Ray->Initial, Lathe->Trans);
  322.  
  323.   MInvTransDirection(D, Ray->Direction, Lathe->Trans);
  324.  
  325.   VLength(len, D);
  326.  
  327.   VInverseScaleEq(D, len);
  328.  
  329.   /* Test if ray misses lathe's cylindrical bound. */
  330.  
  331. #ifdef LATHE_EXTRA_STATS
  332.   Increase_Counter(stats[Lathe_Bound_Tests]);
  333. #endif
  334.  
  335.   if (((D[Y] >= 0.0) && (P[Y] >  Lathe->Height2)) ||
  336.       ((D[Y] <= 0.0) && (P[Y] <  Lathe->Height1)) ||
  337.       ((D[X] >= 0.0) && (P[X] >  Lathe->Radius2)) ||
  338.       ((D[X] <= 0.0) && (P[X] < -Lathe->Radius2)))
  339.   {
  340.     return(FALSE);
  341.   }
  342.  
  343.   /* Get distance r0 of the ray from rotation axis (= y axis). */
  344.  
  345.   r0 = fabs(P[X] * D[Z] - P[Z] * D[X]);
  346.  
  347.   r = D[X] * D[X] + D[Z] * D[Z];
  348.  
  349.   if (r > 0.0)
  350.   {
  351.     r0 /= sqrt(r);
  352.   }
  353.  
  354.   /* Test if ray misses lathe's cylindrical bound. */
  355.  
  356.   if (r0 > Lathe->Radius2)
  357.   {
  358.     return(FALSE);
  359.   }
  360.  
  361. #ifdef LATHE_EXTRA_STATS
  362.   Increase_Counter(stats[Lathe_Bound_Tests_Succeeded]);
  363. #endif
  364.  
  365.   /* Init number of intersections found. */
  366.  
  367.   i = 0;
  368.  
  369.   /* Precalculate some constants that are ray-dependant only. */
  370.  
  371.   m = D[X] * P[X] + D[Z] * P[Z];
  372.  
  373.   Dy2 = D[Y] * D[Y];
  374.  
  375.   /* Intersect all spline segments. */
  376.  
  377.   for (j = 0; j < Lathe->Number; j++)
  378.   {
  379.     Entry = Lathe->Spline->Entry[j];
  380.  
  381.     /* Test if ray misses segment's cylindrical bound. */
  382.  
  383. #ifdef LATHE_EXTRA_STATS
  384.     Increase_Counter(stats[Lathe_Bound_Tests]);
  385. #endif
  386.  
  387.     if (r0 > Entry.r2)
  388.     {
  389.       continue;
  390.     }
  391.  
  392.     if (((D[Y] >= 0.0) && (P[Y] >  Entry.h2)) ||
  393.         ((D[Y] <= 0.0) && (P[Y] <  Entry.h1)) ||
  394.         ((D[X] >= 0.0) && (P[X] >  Entry.r2)) ||
  395.         ((D[X] <= 0.0) && (P[X] < -Entry.r2)))
  396.     {
  397.       continue;
  398.     }
  399.  
  400.     /* Number of roots found. */
  401.  
  402.     n1 = 0;
  403.  
  404.     switch (Lathe->Spline_Type)
  405.     {
  406.       /***********************************************************************
  407.       * Linear spline.
  408.       ************************************************************************/
  409.  
  410.       case LINEAR_SPLINE:
  411.  
  412. #ifdef LATHE_EXTRA_STATS
  413.         Increase_Counter(stats[Lathe_Bound_Tests_Succeeded]);
  414. #endif
  415.  
  416.         /* Solve 2th order polynomial. */
  417.  
  418.         x1[0] = Entry.C[Y] * Entry.C[Y] * r - Entry.C[X] * Entry.C[X] * Dy2;
  419.  
  420.         x1[1] = 2.0 * (Entry.C[Y] * ((Entry.D[Y] - P[Y]) * r + D[Y] * m) - Entry.C[X] * Entry.D[X] * Dy2);
  421.  
  422.         x1[2] = (Entry.D[Y] - P[Y]) * ((Entry.D[Y] - P[Y]) * r + 2.0 * D[Y] * m) + Dy2 * (P[X] * P[X] + P[Z] * P[Z] - Entry.D[X] * Entry.D[X]);
  423.  
  424.         n1 = Solve_Polynomial(2, x1, y1, FALSE, 0.0);
  425.  
  426.         break;
  427.  
  428.       /***********************************************************************
  429.       * Quadratic spline.
  430.       ************************************************************************/
  431.  
  432.       case QUADRATIC_SPLINE:
  433.  
  434.         /* First test the bounding cylinder. */
  435.  
  436.         if (Test_Thick_Cylinder(P, D, Entry.h1, Entry.h2, Sqr(Entry.r1), Sqr(Entry.r2)))
  437.         {
  438. #ifdef LATHE_EXTRA_STATS
  439.           Increase_Counter(stats[Lathe_Bound_Tests_Succeeded]);
  440. #endif
  441.  
  442.           /* Solve 4th order polynomial. */
  443.  
  444.           x1[0] = Entry.B[Y] * Entry.B[Y] * r - Entry.B[X] * Entry.B[X] * Dy2;
  445.  
  446.           x1[1] = 2.0 * (Entry.B[Y] * Entry.C[Y] * r - Entry.B[X] * Entry.C[X] * Dy2);
  447.  
  448.           x1[2] = r * (2.0 * Entry.B[Y] * (Entry.D[Y] - P[Y]) + Entry.C[Y] * Entry.C[Y]) + 2.0 * Entry.B[Y] * D[Y] * m - (2.0 * Entry.B[X] * Entry.D[X] + Entry.C[X] * Entry.C[X]) * Dy2;
  449.  
  450.           x1[3] = 2.0 * (Entry.C[Y] * ((Entry.D[Y] - P[Y]) * r + D[Y] * m) - Entry.C[X] * Entry.D[X] * Dy2);
  451.  
  452.           x1[4] = (Entry.D[Y] - P[Y]) * ((Entry.D[Y] - P[Y]) * r + 2.0 * D[Y] * m) + Dy2 * (P[X] * P[X] + P[Z] * P[Z] - Entry.D[X] * Entry.D[X]);
  453.  
  454.           n1 = Solve_Polynomial(4, x1, y1, Test_Flag(Lathe, STURM_FLAG), 0.0);
  455.         }
  456.  
  457.         break;
  458.  
  459.       /***********************************************************************
  460.       * Cubic spline.
  461.       ************************************************************************/
  462.  
  463.       case CUBIC_SPLINE:
  464.  
  465.         /* First test the bounding cylinder. */
  466.  
  467.         if (Test_Thick_Cylinder(P, D, Entry.h1, Entry.h2, Sqr(Entry.r1), Sqr(Entry.r2)))
  468.         {
  469. #ifdef LATHE_EXTRA_STATS
  470.           Increase_Counter(stats[Lathe_Bound_Tests_Succeeded]);
  471. #endif
  472.  
  473.           /* Solve 6th order polynomial. */
  474.  
  475.           x1[0] = Entry.A[Y] * Entry.A[Y] * r - Entry.A[X] * Entry.A[X] * Dy2;
  476.  
  477.           x1[1] = 2.0 * (Entry.A[Y] * Entry.B[Y] * r - Entry.A[X] * Entry.B[X] * Dy2);
  478.  
  479.           x1[2] = (2.0 * Entry.A[Y] * Entry.C[Y] + Entry.B[Y] * Entry.B[Y]) * r - (2.0 * Entry.A[X] * Entry.C[X] + Entry.B[X] * Entry.B[X]) * Dy2;
  480.  
  481.           x1[3] = 2.0 * ((Entry.A[Y] * Entry.D[Y] + Entry.B[Y] * Entry.C[Y] - Entry.A[Y] * P[Y]) * r + Entry.A[Y] * D[Y] * m - (Entry.A[X] * Entry.D[X] + Entry.B[X] * Entry.C[X]) * Dy2);
  482.  
  483.           x1[4] = (2.0 * Entry.B[Y] * (Entry.D[Y] - P[Y]) + Entry.C[Y] * Entry.C[Y]) * r + 2.0 * Entry.B[Y] * D[Y] * m - (2.0 * Entry.B[X] * Entry.D[X] + Entry.C[X] * Entry.C[X]) * Dy2;
  484.  
  485.           x1[5] = 2.0 * (Entry.C[Y] * ((Entry.D[Y] - P[Y]) * r + D[Y] * m) - Entry.C[X] * Entry.D[X] * Dy2);
  486.  
  487.           x1[6] = (Entry.D[Y] - P[Y]) * ((Entry.D[Y] - P[Y]) * r + 2.0 * D[Y] * m) + Dy2 * (P[X] * P[X] + P[Z] * P[Z] - Entry.D[X] * Entry.D[X]);
  488.  
  489.           n1 = Solve_Polynomial(6, x1, y1, Test_Flag(Lathe, STURM_FLAG), 0.0);
  490.         }
  491.  
  492.         break;
  493.     }
  494.  
  495.     /* Test roots for valid intersections. */
  496.  
  497.     while (n1--)
  498.     {
  499.       w = y1[n1];
  500.  
  501.       if ((w >= 0.0) && (w <= 1.0))
  502.       {
  503.         if (fabs(D[Y]) > EPSILON)
  504.         {
  505.           k = (w * (w * (w * Entry.A[Y] + Entry.B[Y]) + Entry.C[Y]) + Entry.D[Y] - P[Y]) / D[Y];
  506.  
  507.           Intersection[i].n   = j;
  508.           Intersection[i].w   = w;
  509.           Intersection[i++].d = k / len;
  510.  
  511.           if (i > MAX_INTERSECTIONS_PER_SEGMENT * Lathe->Number)
  512.           {
  513.             Error("Too many intersections in intersect_lathe() found.\n");
  514.           }
  515.         }
  516.         else
  517.         {
  518.           k = w * (w * (w * Entry.A[X] + Entry.B[X]) + Entry.C[X]) + Entry.D[X];
  519.  
  520.           x2[0] = r;
  521.           x2[1] = 2.0 * m;
  522.  
  523.           x2[2] = P[X] * P[X] + P[Z] * P[Z] - k * k;
  524.  
  525.           n2 = Solve_Polynomial(2, x2, y2, FALSE, 0.0);
  526.  
  527.           while (n2--)
  528.           {
  529.         k = y2[n2] / len;
  530.  
  531.         /* Test if we already have this intersection. */
  532.  
  533.         for (ii = skip = 0; ii < i; ii++)
  534.         {
  535.           if (fabs(Intersection[ii].d - k) < EPSILON)
  536.           {
  537.             skip = TRUE;
  538.  
  539.         break;
  540.               }
  541.         }
  542.  
  543.         /* Use it if we don't. */
  544.  
  545.         if (!skip)
  546.         {
  547.               Intersection[i].n   = j;
  548.               Intersection[i].w   = w;
  549.               Intersection[i++].d = k;
  550.  
  551.               if (i > MAX_INTERSECTIONS_PER_SEGMENT * Lathe->Number)
  552.               {
  553.                 Error("Too many intersections in intersect_lathe() found.\n");
  554.               }
  555.         }
  556.           }
  557.         }
  558.       }
  559.     }
  560.   }
  561.  
  562.   if (i)
  563.   {
  564.     Increase_Counter(stats[Ray_Lathe_Tests_Succeeded]);
  565.   }
  566.  
  567.   return(i);
  568. }
  569.  
  570.  
  571.  
  572. /*****************************************************************************
  573. *
  574. * FUNCTION
  575. *
  576. *   Inside_Lathe
  577. *
  578. * INPUT
  579. *
  580. *   IPoint - Intersection point
  581. *   Object - Object
  582. *   
  583. * OUTPUT
  584. *   
  585. * RETURNS
  586. *
  587. *   int - TRUE if inside
  588. *   
  589. * AUTHOR
  590. *
  591. *   Dieter Bayer
  592. *   
  593. * DESCRIPTION
  594. *
  595. *   Test if a point lies inside the lathe.
  596. *
  597. * CHANGES
  598. *
  599. *   Jun 1994 : Creation.
  600. *
  601. ******************************************************************************/
  602.  
  603. static int Inside_Lathe(IPoint, Object)
  604. VECTOR IPoint;
  605. OBJECT *Object;
  606. {
  607.   int i, n, NC;
  608.   DBL r, k, w;
  609.   DBL x[4], y[3];
  610.   VECTOR P;
  611.   LATHE_SPLINE_ENTRY Entry;
  612.   LATHE *Lathe = (LATHE *)Object;
  613.  
  614.   /* Transform the point into the lathe space. */
  615.  
  616.   MInvTransPoint(P, IPoint, Lathe->Trans);
  617.  
  618.   /* Number of crossings. */
  619.  
  620.   NC = 0;
  621.  
  622.   if ((P[Y] >= Lathe->Height1) && (P[Y] <= Lathe->Height2))
  623.   {
  624.     r = sqrt(P[X] * P[X] + P[Z] * P[Z]);
  625.  
  626.     if ((r >= Lathe->Radius1) && (r <= Lathe->Radius2))
  627.     {
  628.       for (i = 0; i < Lathe->Number; i++)
  629.       {
  630.         Entry = Lathe->Spline->Entry[i];
  631.  
  632.         /* Test if we are inside the segments cylindrical bound. */
  633.  
  634.         if ((P[Y] - Entry.h1 >= -EPSILON) &&
  635.             (Entry.h2 - P[Y] >= -EPSILON) && (r <= Entry.r2))
  636.         {
  637.           x[0] = Entry.A[Y];
  638.           x[1] = Entry.B[Y];
  639.           x[2] = Entry.C[Y];
  640.           x[3] = Entry.D[Y] - P[Y];
  641.  
  642.           n = Solve_Polynomial(3, x, y, Test_Flag(Lathe, STURM_FLAG), 0.0);
  643.  
  644.           while (n--)
  645.           {
  646.             w = y[n];
  647.  
  648.             if ((w >= 0.0) && (w <= 1.0))
  649.             {
  650.               k  = w * (w * (w * Entry.A[X] + Entry.B[X]) + Entry.C[X]) + Entry.D[X] - r;
  651.  
  652.               if (k >= 0.0)
  653.               {
  654.                 NC++;
  655.               }
  656.             }
  657.           }
  658.         }
  659.       }
  660.     }
  661.   }
  662.  
  663.   if (NC & 1)
  664.   {
  665.     return(!Test_Flag(Lathe, INVERTED_FLAG));
  666.   }
  667.   else
  668.   {
  669.     return(Test_Flag(Lathe, INVERTED_FLAG));
  670.   }
  671. }
  672.  
  673.  
  674.  
  675. /*****************************************************************************
  676. *
  677. * FUNCTION
  678. *
  679. *   Lathe_Normal
  680. *
  681. * INPUT
  682. *
  683. *   Result - Normal vector
  684. *   Object - Object
  685. *   Inter  - Intersection found
  686. *   
  687. * OUTPUT
  688. *
  689. *   Result
  690. *   
  691. * RETURNS
  692. *   
  693. * AUTHOR
  694. *
  695. *   Dieter Bayer
  696. *   
  697. * DESCRIPTION
  698. *
  699. *   Calculate the normal of the lathe in a given point.
  700. *
  701. * CHANGES
  702. *
  703. *   Jun 1994 : Creation.
  704. *
  705. ******************************************************************************/
  706.  
  707. static void Lathe_Normal(Result, Object, Inter)
  708. OBJECT *Object;
  709. VECTOR Result;
  710. INTERSECTION *Inter;
  711. {
  712.   DBL r, dx, dy;
  713.   VECTOR P, N;
  714.   LATHE *Lathe = (LATHE *)Object;
  715.   LATHE_SPLINE_ENTRY Entry;
  716.  
  717.   Entry = Lathe->Spline->Entry[Inter->i1];
  718.  
  719.   /* Transform the point into the lathe space. */
  720.  
  721.   MInvTransPoint(P, Inter->IPoint, Lathe->Trans);
  722.  
  723.   /* Get distance from rotation axis. */
  724.  
  725.   r = P[X] * P[X] + P[Z] * P[Z];
  726.  
  727.   if (r > EPSILON)
  728.   {
  729.     r = sqrt(r);
  730.  
  731.     /* Get derivatives. */
  732.  
  733.     dx = Inter->d1 * (3.0 * Entry.A[X] * Inter->d1 + 2.0 * Entry.B[X]) + Entry.C[X];
  734.     dy = Inter->d1 * (3.0 * Entry.A[Y] * Inter->d1 + 2.0 * Entry.B[Y]) + Entry.C[Y];
  735.  
  736.     /* Get normal by rotation. */
  737.  
  738.     N[X] =  dy * P[X];
  739.     N[Y] = -dx * r;
  740.     N[Z] =  dy * P[Z];
  741.   }
  742.   else
  743.   {
  744.     N[X] = N[Z] = 0.0; N[Y] = 1.0;
  745.   }
  746.  
  747.   /* Transform the normalt out of the lathe space. */
  748.  
  749.   MTransNormal(Result, N, Lathe->Trans);
  750.  
  751.   VNormalize(Result, Result);
  752. }
  753.  
  754.  
  755.  
  756. /*****************************************************************************
  757. *
  758. * FUNCTION
  759. *
  760. *   Translate_Lathe
  761. *
  762. * INPUT
  763. *
  764. *   Object - Object
  765. *   Vector - Translation vector
  766. *   
  767. * OUTPUT
  768. *
  769. *   Object
  770. *   
  771. * RETURNS
  772. *   
  773. * AUTHOR
  774. *
  775. *   Dieter Bayer
  776. *   
  777. * DESCRIPTION
  778. *
  779. *   Translate a lathe.
  780. *
  781. * CHANGES
  782. *
  783. *   Jun 1994 : Creation.
  784. *
  785. ******************************************************************************/
  786.  
  787. static void Translate_Lathe(Object, Vector, Trans)
  788. OBJECT *Object;
  789. VECTOR Vector;
  790. TRANSFORM *Trans;
  791. {
  792.   Transform_Lathe(Object, Trans);
  793. }
  794.  
  795.  
  796.  
  797. /*****************************************************************************
  798. *
  799. * FUNCTION
  800. *
  801. *   Rotate_Lathe
  802. *
  803. * INPUT
  804. *
  805. *   Object - Object
  806. *   Vector - Rotation vector
  807. *   
  808. * OUTPUT
  809. *
  810. *   Object
  811. *   
  812. * RETURNS
  813. *   
  814. * AUTHOR
  815. *
  816. *   Dieter Bayer
  817. *   
  818. * DESCRIPTION
  819. *
  820. *   Rotate a lathe.
  821. *
  822. * CHANGES
  823. *
  824. *   Jun 1994 : Creation.
  825. *
  826. ******************************************************************************/
  827.  
  828. static void Rotate_Lathe(Object, Vector, Trans)
  829. OBJECT *Object;
  830. VECTOR Vector;
  831. TRANSFORM *Trans;
  832. {
  833.   Transform_Lathe(Object, Trans);
  834. }
  835.  
  836.  
  837.  
  838. /*****************************************************************************
  839. *
  840. * FUNCTION
  841. *
  842. *   Scale_Lathe
  843. *
  844. * INPUT
  845. *
  846. *   Object - Object
  847. *   Vector - Scaling vector
  848. *   
  849. * OUTPUT
  850. *
  851. *   Object
  852. *   
  853. * RETURNS
  854. *   
  855. * AUTHOR
  856. *
  857. *   Dieter Bayer
  858. *   
  859. * DESCRIPTION
  860. *
  861. *   Scale a lathe.
  862. *
  863. * CHANGES
  864. *
  865. *   Jun 1994 : Creation.
  866. *
  867. ******************************************************************************/
  868.  
  869. static void Scale_Lathe(Object, Vector, Trans)
  870. OBJECT *Object;
  871. VECTOR Vector;
  872. TRANSFORM *Trans;
  873. {
  874.   Transform_Lathe(Object, Trans);
  875. }
  876.  
  877.  
  878.  
  879. /*****************************************************************************
  880. *
  881. * FUNCTION
  882. *
  883. *   Transform_Lathe
  884. *
  885. * INPUT
  886. *
  887. *   Object - Object
  888. *   Trans  - Transformation to apply
  889. *   
  890. * OUTPUT
  891. *
  892. *   Object
  893. *   
  894. * RETURNS
  895. *   
  896. * AUTHOR
  897. *
  898. *   Dieter Bayer
  899. *   
  900. * DESCRIPTION
  901. *
  902. *   Transform a lathe and recalculate its bounding box.
  903. *
  904. * CHANGES
  905. *
  906. *   Jun 1994 : Creation.
  907. *
  908. ******************************************************************************/
  909.  
  910. static void Transform_Lathe(Object, Trans)
  911. OBJECT *Object;
  912. TRANSFORM *Trans;
  913. {
  914.   Compose_Transforms(((LATHE *)Object)->Trans, Trans);
  915.  
  916.   Compute_Lathe_BBox((LATHE *)Object);
  917. }
  918.  
  919.  
  920.  
  921. /*****************************************************************************
  922. *
  923. * FUNCTION
  924. *
  925. *   Invert_Lathe
  926. *
  927. * INPUT
  928. *
  929. *   Object - Object
  930. *   
  931. * OUTPUT
  932. *
  933. *   Object
  934. *   
  935. * RETURNS
  936. *   
  937. * AUTHOR
  938. *
  939. *   Dieter Bayer
  940. *   
  941. * DESCRIPTION
  942. *
  943. *   Invert a lathe.
  944. *
  945. * CHANGES
  946. *
  947. *   Jun 1994 : Creation.
  948. *
  949. ******************************************************************************/
  950.  
  951. static void Invert_Lathe(Object)
  952. OBJECT *Object;
  953. {
  954.   Invert_Flag(Object, INVERTED_FLAG);
  955. }
  956.  
  957.  
  958.  
  959. /*****************************************************************************
  960. *
  961. * FUNCTION
  962. *
  963. *   Create_Lathe
  964. *
  965. * INPUT
  966. *   
  967. * OUTPUT
  968. *   
  969. * RETURNS
  970. *
  971. *   LATHE * - new lathe
  972. *   
  973. * AUTHOR
  974. *
  975. *   Dieter Bayer
  976. *   
  977. * DESCRIPTION
  978. *
  979. *   Create a new lathe.
  980. *
  981. * CHANGES
  982. *
  983. *   Jun 1994 : Creation.
  984. *
  985. ******************************************************************************/
  986.  
  987. LATHE *Create_Lathe()
  988. {
  989.   LATHE *New;
  990.  
  991.   New = (LATHE *)POV_MALLOC(sizeof(LATHE), "lathe");
  992.  
  993.   INIT_OBJECT_FIELDS(New,LATHE_OBJECT,&Lathe_Methods)
  994.  
  995.   New->Trans = Create_Transform();
  996.  
  997.   New->Spline_Type = LINEAR_SPLINE;
  998.  
  999.   New->Number = 0;
  1000.  
  1001.   New->Spline = NULL;
  1002.  
  1003.   New->Radius1 =
  1004.   New->Radius2 =
  1005.   New->Height1 =
  1006.   New->Height2 = 0.0;
  1007.  
  1008.   return(New);
  1009. }
  1010.  
  1011.  
  1012.  
  1013. /*****************************************************************************
  1014. *
  1015. * FUNCTION
  1016. *
  1017. *   Copy_Lathe
  1018. *
  1019. * INPUT
  1020. *
  1021. *   Object - Object
  1022. *   
  1023. * OUTPUT
  1024. *   
  1025. * RETURNS
  1026. *
  1027. *   void * - New lathe
  1028. *   
  1029. * AUTHOR
  1030. *
  1031. *   Dieter Bayer
  1032. *   
  1033. * DESCRIPTION
  1034. *
  1035. *   Copy a lathe structure.
  1036. *
  1037. *   NOTE: The splines are not copied, only the number of references is
  1038. *         counted, so that Destray_Lathe() knows if they can be destroyed.
  1039. *
  1040. * CHANGES
  1041. *
  1042. *   Jun 1994 : Creation.
  1043. *
  1044. *   Sep 1994 : Fixed memory leakage bug. [DB]
  1045. *
  1046. ******************************************************************************/
  1047.  
  1048. static void *Copy_Lathe(Object)
  1049. OBJECT *Object;
  1050. {
  1051.   LATHE *New, *Lathe = (LATHE *)Object;
  1052.  
  1053.   New = Create_Lathe();
  1054.  
  1055.   /* Get rid of the transformation created in Create_Lathe(). */
  1056.  
  1057.   Destroy_Transform(New->Trans);
  1058.  
  1059.   /* Copy lathe. */
  1060.  
  1061.   *New = *Lathe;
  1062.  
  1063.   New->Trans = Copy_Transform(Lathe->Trans);
  1064.  
  1065.   New->Spline->References++;
  1066.  
  1067.   return(New);
  1068. }
  1069.  
  1070.  
  1071.  
  1072. /*****************************************************************************
  1073. *
  1074. * FUNCTION
  1075. *
  1076. *   Destroy_Lathe
  1077. *
  1078. * INPUT
  1079. *
  1080. *   Object - Object
  1081. *   
  1082. * OUTPUT
  1083. *
  1084. *   Object
  1085. *   
  1086. * RETURNS
  1087. *   
  1088. * AUTHOR
  1089. *
  1090. *   Dieter Bayer
  1091. *   
  1092. * DESCRIPTION
  1093. *
  1094. *   Destroy a lathe.
  1095. *
  1096. *   NOTE: The splines are destroyed if they are no longer used by any copy.
  1097. *
  1098. * CHANGES
  1099. *
  1100. *   Jun 1994 : Creation.
  1101. *
  1102. ******************************************************************************/
  1103.  
  1104. static void Destroy_Lathe(Object)
  1105. OBJECT *Object;
  1106. {
  1107.   LATHE *Lathe = (LATHE *)Object;
  1108.  
  1109.   Destroy_Transform(Lathe->Trans);
  1110.  
  1111.   if (--(Lathe->Spline->References) == 0)
  1112.   {
  1113.     POV_FREE (Lathe->Spline->Entry);
  1114.     POV_FREE (Lathe->Spline);
  1115.   }
  1116.  
  1117.   POV_FREE (Object);
  1118. }
  1119.  
  1120.  
  1121.  
  1122. /*****************************************************************************
  1123. *
  1124. * FUNCTION
  1125. *
  1126. *   Compute_Lathe_BBox
  1127. *
  1128. * INPUT
  1129. *
  1130. *   Lathe - Lathe
  1131. *   
  1132. * OUTPUT
  1133. *
  1134. *   Lathe
  1135. *   
  1136. * RETURNS
  1137. *   
  1138. * AUTHOR
  1139. *
  1140. *   Dieter Bayer
  1141. *   
  1142. * DESCRIPTION
  1143. *
  1144. *   Calculate the bounding box of a lathe.
  1145. *
  1146. * CHANGES
  1147. *
  1148. *   Jun 1994 : Creation.
  1149. *
  1150. ******************************************************************************/
  1151.  
  1152. void Compute_Lathe_BBox(Lathe)
  1153. LATHE *Lathe;
  1154. {
  1155.   Make_BBox(Lathe->BBox, -Lathe->Radius2, Lathe->Height1, -Lathe->Radius2,
  1156.     2.0 * Lathe->Radius2, Lathe->Height2 - Lathe->Height1, 2.0 * Lathe->Radius2);
  1157.  
  1158.   Recompute_BBox(&Lathe->BBox, Lathe->Trans);
  1159. }
  1160.  
  1161.  
  1162.  
  1163. /*****************************************************************************
  1164. *
  1165. * FUNCTION
  1166. *
  1167. *   Compute_Lathe
  1168. *
  1169. * INPUT
  1170. *
  1171. *   Lathe - Lathe
  1172. *   P     - Points defining lathe
  1173. *   
  1174. * OUTPUT
  1175. *
  1176. *   Lathe
  1177. *   
  1178. * RETURNS
  1179. *   
  1180. * AUTHOR
  1181. *
  1182. *   Dieter Bayer
  1183. *   
  1184. * DESCRIPTION
  1185. *
  1186. *   Calculate the spline segments of a lathe from a set of points.
  1187. *
  1188. *   Note that the number of points in the lathe has to be set.
  1189. *
  1190. * CHANGES
  1191. *
  1192. *   Jun 1994 : Creation.
  1193. *
  1194. ******************************************************************************/
  1195.  
  1196. void Compute_Lathe(Lathe, P)
  1197. LATHE *Lathe;
  1198. UV_VECT *P;
  1199. {
  1200.   int i, n;
  1201.   DBL x[4], xmin, xmax;
  1202.   DBL y[4], ymin, ymax;
  1203.   DBL c[3], r[2];
  1204.   UV_VECT A, B, C, D;
  1205.  
  1206.   /* Allocate Lathe->Number segments. */
  1207.  
  1208.   if (Lathe->Spline == NULL)
  1209.   {
  1210.     Lathe->Spline = (LATHE_SPLINE *)POV_MALLOC(sizeof(LATHE_SPLINE), "spline segments of lathe");
  1211.  
  1212.     /* Init spline. */
  1213.  
  1214.     Lathe->Spline->References = 1;
  1215.  
  1216.     Lathe->Spline->Entry = (LATHE_SPLINE_ENTRY *)POV_MALLOC(Lathe->Number*sizeof(LATHE_SPLINE_ENTRY), "spline segments of lathe");
  1217.   }
  1218.   else
  1219.   {
  1220.     /* This should never happen! */
  1221.  
  1222.     Error("Lathe segments are already defined.\n");
  1223.   }
  1224.  
  1225.   /***************************************************************************
  1226.   * Calculate segments.
  1227.   ****************************************************************************/
  1228.  
  1229.   /* We want to know the size of the overall bounding cylinder. */
  1230.  
  1231.   xmax = ymax = -BOUND_HUGE;
  1232.   xmin = ymin = BOUND_HUGE;
  1233.  
  1234.   for (i = 0; i < Lathe->Number; i++)
  1235.   {
  1236.     switch (Lathe->Spline_Type)
  1237.     {
  1238.       /*************************************************************************
  1239.       * Linear spline (nothing more than a simple polygon).
  1240.       **************************************************************************/
  1241.  
  1242.       case LINEAR_SPLINE:
  1243.  
  1244.         /* Use linear interpolation. */
  1245.  
  1246.         A[X] = 0.0;
  1247.         B[X] = 0.0;
  1248.         C[X] = -1.0 * P[i][X] + 1.0 * P[i+1][X];
  1249.         D[X] =  1.0 * P[i][X];
  1250.  
  1251.         A[Y] = 0.0;
  1252.         B[Y] = 0.0;
  1253.         C[Y] = -1.0 * P[i][Y] + 1.0 * P[i+1][Y];
  1254.         D[Y] =  1.0 * P[i][Y];
  1255.  
  1256.         /* Get maximum coordinates in current segment. */
  1257.  
  1258.         x[0] = x[2] = P[i][X];
  1259.         x[1] = x[3] = P[i+1][X];
  1260.  
  1261.         y[0] = y[2] = P[i][Y];
  1262.         y[1] = y[3] = P[i+1][Y];
  1263.  
  1264.         break;
  1265.  
  1266.  
  1267.       /*************************************************************************
  1268.       * Quadratic spline.
  1269.       **************************************************************************/
  1270.  
  1271.       case QUADRATIC_SPLINE:
  1272.  
  1273.         /* Use quadratic interpolation. */
  1274.  
  1275.         A[X] =  0.0;
  1276.         B[X] =  0.5 * P[i][X] - 1.0 * P[i+1][X] + 0.5 * P[i+2][X];
  1277.         C[X] = -0.5 * P[i][X]                   + 0.5 * P[i+2][X];
  1278.         D[X] =                  1.0 * P[i+1][X];
  1279.  
  1280.         A[Y] =  0.0;
  1281.         B[Y] =  0.5 * P[i][Y] - 1.0 * P[i+1][Y] + 0.5 * P[i+2][Y];
  1282.         C[Y] = -0.5 * P[i][Y]                   + 0.5 * P[i+2][Y];
  1283.         D[Y] =                  1.0 * P[i+1][Y];
  1284.  
  1285.         /* Get maximum coordinates in current segment. */
  1286.  
  1287.         x[0] = x[2] = P[i+1][X];
  1288.         x[1] = x[3] = P[i+2][X];
  1289.  
  1290.         y[0] = y[2] = P[i+1][Y];
  1291.         y[1] = y[3] = P[i+2][Y];
  1292.  
  1293.         break;
  1294.  
  1295.  
  1296.       /*************************************************************************
  1297.       * Cubic spline.
  1298.       **************************************************************************/
  1299.  
  1300.       case CUBIC_SPLINE:
  1301.  
  1302.         /* Use cubic interpolation. */
  1303.  
  1304.         A[X] = -0.5 * P[i][X] + 1.5 * P[i+1][X] - 1.5 * P[i+2][X] + 0.5 * P[i+3][X];
  1305.         B[X] =        P[i][X] - 2.5 * P[i+1][X] + 2.0 * P[i+2][X] - 0.5 * P[i+3][X];
  1306.         C[X] = -0.5 * P[i][X]                   + 0.5 * P[i+2][X];
  1307.         D[X] =                        P[i+1][X];
  1308.  
  1309.         A[Y] = -0.5 * P[i][Y] + 1.5 * P[i+1][Y] - 1.5 * P[i+2][Y] + 0.5 * P[i+3][Y];
  1310.         B[Y] =        P[i][Y] - 2.5 * P[i+1][Y] + 2.0 * P[i+2][Y] - 0.5 * P[i+3][Y];
  1311.         C[Y] = -0.5 * P[i][Y]                   + 0.5 * P[i+2][Y];
  1312.         D[Y] =                        P[i+1][Y];
  1313.  
  1314.         /* Get maximum coordinates in current segment. */
  1315.  
  1316.         x[0] = x[2] = P[i+1][X];
  1317.         x[1] = x[3] = P[i+2][X];
  1318.  
  1319.         y[0] = y[2] = P[i+1][Y];
  1320.         y[1] = y[3] = P[i+2][Y];
  1321.  
  1322.         break;
  1323.  
  1324.       default:
  1325.  
  1326.         Error("Unknown lathe type in Compute_Lathe().\n");
  1327.  
  1328.     }
  1329.  
  1330.     Assign_UV_Vect(Lathe->Spline->Entry[i].A, A);
  1331.     Assign_UV_Vect(Lathe->Spline->Entry[i].B, B);
  1332.     Assign_UV_Vect(Lathe->Spline->Entry[i].C, C);
  1333.     Assign_UV_Vect(Lathe->Spline->Entry[i].D, D);
  1334.  
  1335.     /* Get maximum coordinates in current segment. */
  1336.  
  1337.     c[0] = 3.0 * A[X];
  1338.     c[1] = 2.0 * B[X];
  1339.     c[2] = C[X];
  1340.  
  1341.     n = Solve_Polynomial(2, c, r, FALSE, 0.0);
  1342.  
  1343.     while (n--)
  1344.     {
  1345.       if ((r[n] >= 0.0) && (r[n] <= 1.0))
  1346.       {
  1347.         x[n] = r[n] * (r[n] * (r[n] * A[X] + B[X]) + C[X]) + D[X];
  1348.       }
  1349.     }
  1350.  
  1351.     c[0] = 3.0 * A[Y];
  1352.     c[1] = 2.0 * B[Y];
  1353.     c[2] = C[Y];
  1354.  
  1355.     n = Solve_Polynomial(2, c, r, FALSE, 0.0);
  1356.  
  1357.     while (n--)
  1358.     {
  1359.       if ((r[n] >= 0.0) && (r[n] <= 1.0))
  1360.       {
  1361.         y[n] = r[n] * (r[n] * (r[n] * A[Y] + B[Y]) + C[Y]) + D[Y];
  1362.       }
  1363.     }
  1364.  
  1365.     /* Set current segment's bounding cylinder. */
  1366.  
  1367.     Lathe->Spline->Entry[i].r1 = min(min(x[0], x[1]), min(x[2], x[3]));
  1368.     Lathe->Spline->Entry[i].r2 = max(max(x[0], x[1]), max(x[2], x[3]));
  1369.  
  1370.     Lathe->Spline->Entry[i].h1 = min(min(y[0], y[1]), min(y[2], y[3]));
  1371.     Lathe->Spline->Entry[i].h2 = max(max(y[0], y[1]), max(y[2], y[3]));
  1372.  
  1373.     /* Keep track of overall bounding cylinder. */
  1374.  
  1375.     xmin = min(xmin, Lathe->Spline->Entry[i].r1);
  1376.     xmax = max(xmax, Lathe->Spline->Entry[i].r2);
  1377.  
  1378.     ymin = min(ymin, Lathe->Spline->Entry[i].h1);
  1379.     ymax = max(ymax, Lathe->Spline->Entry[i].h2);
  1380.   }
  1381.  
  1382.   /* Set overall bounding cylinder. */
  1383.  
  1384.   Lathe->Radius1 = xmin;
  1385.   Lathe->Radius2 = xmax;
  1386.  
  1387.   Lathe->Height1 = ymin;
  1388.   Lathe->Height2 = ymax;
  1389. }
  1390.  
  1391.  
  1392.  
  1393.